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ABSTRACT 

We present the results of axisymmetric time-dependent hydrodynamic calculations of 
line-driven winds from accretion disks in active galactic nuclei (AGN). We assume the disk 
is fiat, Keplerian, geometrically thin, and optically thick, radiating according to the a-disk 
prescription. The central engine of the AGN is a source of both ionizing X-rays and wind-driving 
ultraviolet (UV) photons. To calculate the radiation force, we take into account radiation from 
the disk and the central engine. The gas temperature and ionization state in the wind are 
calculated self-consistently from the photoionization and heating rate of the central engine. 

We find that a disk accreting onto a 10 8 M Q black hole at the rate of 1.8 M Q yr _1 can 
launch a wind at ~ 10 16 cm from the central engine. The X-rays from the central object are 
significantly attenuated by the disk atmosphere so they cannot prevent the local disk radiation 
from pushing matter away from the disk. However in the supersonic portion of the flow high 
above the disk, the X-rays can overionize the gas and decrease the wind terminal velocity. For a 
reasonable X-ray opacity, e.g., kx = 40 g _1 cm 2 , the disk wind can be accelerated by the central 
UV radiation to velocities of up to 15000 km s _1 at a distance of ~ 10 17 cm from the central 
engine. The covering factor of the disk wind is ~ 0.2. The wind is unsteady and consists of an 
opaque, slow vertical flow near the disk that is bounded on the polar side by a high-velocity 
stream. A typical column density through the fast stream is a few 10 23 cm -2 so the stream 
is optically thin to the UV radiation. This low column density is precisely why gas can be 
accelerated to high velocities. The fast stream contributes nearly 100% to the total wind mass 
loss rate of 0.5 M yr -1 . 



Subject headings: accretion disks - hydrodynamics - instabilities - outflows - galaxies: active - 
methods: numerical 
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Many observed spectral features of active galactic nuclei (AGN) indicate that outflows are common in 
these systems. The prominent broad emission lines (BEL) in ultraviolet (UV) from H I, O VI, N V, C IV, 
and Si IV are the defining feature of quasars (Blandford, Netzer & Woltjer 1990; Osterbrock 1989), and they 
may be associated with a high velocity outflow. Fast outflows can also explain narrow UV absorption lines 
from highly ionized species - such as C IV and N V - observed in approximately half of the HST-observed 
Seyfert 1 galaxies (e.g., Crenshaw 1997). These narrow absorption lines are all blucshifted relative to the 
systemic velocity by to -1600 km s _1 . Recent ASCA observations show that in a small sample observed 
by HST and ASCA, all of the Seyfert galaxies with warm absorbers also show intrinsic UV absorption. 
A high-resolution X-ray observation of the Seyfert galaxy NGC 5548 obtained by Chandra shows strong, 
narrow absorption lines from highly ionized species (Kaastra et al. 2000). The lines are blucshifted by a 
few hundred km s -1 and they are reminiscent of the narrow absorption lines observed in the UV. Some 
QSOs also have intrinsic UV narrow-line absorbers with line-of-sight velocities as large as 51000 km s _1 
(Hamann ct al. 1997). Perhaps the most compelling evidence for fast outflows in QSOs are strong broad 
absorption lines (BALs) in the UV resonance lines of highly ionized species such as N V, C IV, and Si IV. 
The BALs are always blueshifted relative to the emission-line rest frame, indicating the presence of outflows 
from the active nucleus, with velocities as large as 60000 kms -1 . 

A key constraint on any model for the origin of AGN outflows is the ionization balance. On one hand 
we observe very high luminosities in X-rays and the UV and on the other hand we observe spectral lines 
from moderately and highly ionized species. One wonders then how the gas avoids full photoionization and 
we see any spectral lines at all. Generally, two mechanisms have been proposed to resolve the so-called 
overionizaton problem: (i) the AGN outflows have filling factors less than one and consist of dense clouds 
and (ii) the filling factors equal one but the outflows are shielded from the powerful radiation by some 
material located between the central engine and the outflow (e.g., Krolik 1999, and below). 

Many theoretical models have been proposed to explain outflows in AGNs. Comprehensive reviews of 
recent theoretical work on outflows in AGNs can be found in Arav, Shlosman & Weymann (1997). Several 
forces have been suggested to accelerate outflows in AGN, for example, gas pressure (e.g., Weymann ct 
al. 1982; Begelman, de Kool & Sikora 1991), magneto-centrifugal force due to an accretion disk (e.g., 
Blandford & Payne 1982; Emmering, Blandford & Shlosman 1992; Konigl & Kartje 1994; Bottorff et al. 
1997), radiation pressure due to dust (e.g., Voit, Weymann & Korista 1993; Scoville & Norman 1995), and 
radiation pressure due to lines (e.g., Drew & Bokscnberg 1984; Shlosman, Vitello & Shaviv 1985; Arav, Li 
& Begelman 1994; dc Kool & Begelman 1995, hereafter dKB; Murray et al. 1995, hereafter MCGV). Most 
of these models avoid overionization of the wind by assuming high density clouds. 

One plausible scenario for the AGN outflows is that radiation pressure on spectral lines drives a wind 
from an accretion disk around a black hole. The presence of the BALs themselves strongly indicates that 
substantial momentum is transfered from a powerful radiation field to the gas. A crucial clue to the origin 
of AGN outflows comes from the discovery of line-locking in BAL QSO spectra (e.g., Foltz et al. 1987). 
Weymann et al. (1991) discovered that a composite spectrum of their BAL QSO sample shows a double 
trough in the C IV A1549 BAL, separated in velocity space by the N V A 1240-Lya splitting ~ 5, 900 km s _1 
(see also Korista et al. 1993). In fact Arav & Begelman (1994, see also Arav et al. 1995 and Arav 1996) 
showed that an absorption hump in the C IV A BAL is "the ghost of Lya" due to modulation of the 
radiation force by the strong emission in Lya. 

Our understanding of how line-driving produces powerful high velocity winds is based on the studies of 
winds in hot stars that radiate mostly in the UV (e.g., Castor, Abbott, Klein 1975, hereafter CAK; Abbott 
1982). The essential concept underpinning these models is that the momentum is extracted most efficiently 
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from the radiation field via line opacity. With the inclusion of lines, CAK showed that the effective radiation 
force can increase by several orders of magnitude above that due to electron-scattering alone. Thus even a 
star that radiates at around 0.1% of its Eddington limit, Lem can have a strong wind. Up to now, however, 
it has been difficult to apply line-driven stellar wind models to QSOs because of the fundamental difference 
in geometry: stellar winds are to a good approximation spherically symmetric, whereas the wind in BAL 
QSOs likely arises from a flattened disk and is therefore axially symmetric. 

Line-radiation driven disk wind models for BAL QSOs have been proposed by several authors. For 
example, MCGV studied a wind arising from a small distance from the central engine (~ 10 16 cm). Their 
model relies on the local disk radiation to launch gas from the disk and on the central radiation to accelerate 
the gas in the radial direction to high velocities. A key ingredient to this model is that the central engine 
radiation is attenuated by 'hitchhiking' gas located between the central engine and wind. Rather than 
showing from first principles the origin of this hitchhiking gas, they simply give plausible arguments for its 
existence. MCGV adopted assumptions that require the flow to be time-independent and restrict the flow 
geometry to a quasi 1-D radial flow. Other models of radiation-driven outflows in AGNs propose that the 
radiation force accelerates the wind that has been launched and kept at a low ionization state by different 
mechanisms. For example, the dKB model supposes that the wind is seated at a large distance (~ 10 18 cm) 
from the central engine. Because the disk is cool at these radii, the gas is not lifted from the disk by 
radiation pressure but by some alternative mechanism. Once the gas is high enough above the disk, the 
radiation force due to the central engine accelerates the gas to high velocity. The gas is not over-ionized in 
this model due to a small filling factor caused by strong magnetic fields which confine the dense clouds. 

Recently it has been possible to model line-driven disk winds using 2-dimensional axisymmetric 
numerical hydrodynamical simulations (e.g., Pereyra, Kallman & Blondin 1997; PSD I; Proga 1999; Proga, 
Stone & Drew 1999, hereafter PSD II). All these models have been calculated for white dwarf accretion 
disks. PSD I found that line-driven disk winds are produced only when the effective luminosity of the 
disk (i.e., the luminosity of the disk times the maximum value of the force multiplier, LvM max ) exceeds 
the Eddington limit, LEdd- If the dominant contribution to the total radiation field comes from the 
disk, then the outflow is intrinsically unsteady and characterized by large amplitude velocity and density 
fluctuations. On the other hand, if the total luminosity of the system is dominated by the central object, 
then the outflow is steady. In either case, PSD I and PSD II found that the structure of the wind consists 
of a dense, slow outflow, that is bounded on the polar side by a high-velocity, lower density stream. The 
flow geometry is controlled largely by the geometry of the radiation field - a brighter disk/central object 
produces a more polar/equatorial wind. Global properties such as the total mass loss rate and terminal 
velocity depend more on the system luminosity and are insensitive to geometry. The mass loss rate is a 
strong function of the effective Eddington luminosity and is of the same order of magnitude as that of a 
simple spherically-symmetric stellar wind (see also Proga 1999). Matter is fed into the fast stream from 
within a few central object radii. The terminal velocity of the stream is similar to that of the terminal 
velocity of a corresponding spherical stellar wind, i.e., ~ a few v esc , where v esc , is the escape velocity 
from the photosphere. 

In this paper, we use the methods developed by PSD to study disk winds in AGNs. However to study 
AGNs winds, we need to relax some of the assumptions adopted by PSD. In particular, we can no longer 
assume the outflowing gas is isothermal, has a fixed ionization state, or is optically thin. Instead, we adopt 
a simplified treatment of photoionization, and radiative cooling and heating that allow us to compute 
self-consistently the ionization state, and therefore line force, in the wind. 

Here we calculate a few disk wind models using our extensions to the PSD II method. We concentrate 
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on assessing how winds can be driven from a disk in the presence of very strong ionizing radiation. We 
describe our calculation in Section 2. We present our results in Section 3 and discuss them together with 
perceived limitations in Section 4. The paper ends, in section 5, with our conclusions. 



2. Method 



2.1. Hydrodynamics 



To calculate the structure and evolution of a wind from a disk, we solve the equations of hydrodynamics 

g + pV.v = 0, (1) 

p^ = -VP + pg + P F rad (2) 
(z) =-pV-v + p£, (3) 



Dt VP, 

where p is the mass density, P is the gas pressure, v is the velocity, e is the internal energy density, C is the 
net cooling rate, g is the gravitational acceleration of the central object, F rad is the total radiation force 
per unit mass. We adopt an adiabatic equation of state P = (7 — l)e, and consider models with 7 = 5/3. 

Our calculations are performed in spherical polar coordinates (r,6,<f)). We assume axial symmetry 
about the rotational axis of the accretion disk (9 = 0°). However the 9 = 90° axis is not coincident with 
the disk midplane. We allow this axis to be above the midplane, at the height, z , corresponding to the disk 
pressure scale height. Figure 1 shows a schematic representation of our computation domain. The offset 
enters our calculations in the computation of the gravity; for example the gravity has non-zero 9 component 
at 6 = 90°, and in the computation of the Keplerian velocity along the 9 = 90° axis. 

This offset is introduced so that the bottom of our computational grid is located at the disk photosphere, 
where the wind is launched. Unlike the thin disk studied by PSD, the radiation dominated disk in AGN 
may have significant geometrical width. If the base of the computational grid is placed at the equatorial 
plane (midplane), this would require modeling the full internal structure of the disk, which is complex and 
turbulent (e.g., Balbus & Hawley 1998). Instead, we locate the base of the grid where the vertical radiation 
force due to electron scattering cancels out the vertical component of gravity. 

A realistic description of the radiation field from the central engine would require detailed knowledge 
of the geometry of the flow near the central engine, which is beyond the scope of this investigation. Instead 
we model the central engine as a point source of radiation located at the origin of our computational grid 
(that is located a height z a above the disk midplane). This implies the central engine has finite width, 
perhaps associated with a hot corona. We expect the high column density of the disk to attenuate radiation 
close to the 9 = 90° plane. 

Our standard computational domain is defined to occupy the radial range = 100 r* < r < 
r D = 1000 r», where r* is the inner radius of the disk, and the angular range 0° < 9 < 90°. For comparison, 
we also calculate some models with the radial range 50 r* < r < 1000 r* and with the radial range 
200 r* < r < 1000 r„. The r — 9 domain is discrctized into zones. Our numerical resolution consists of 100 
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zones in each of the r and 9 directions, with fixed zone size ratios, drk+i/drk — 1-05 and d6i/d6i+i = 1.087. 
Gridding in this manner ensures good spatial resolution close to the radiating surface of the disk and the 
inner boundary at rj. 

The boundary conditions are specified as follows. At 6 = 0, we apply an axis-of-symmetry boundary 
condition. For the outer radial boundary, we apply an outflow boundary condition. For the inner radial 
boundary r — r* and for 9 — 90°, we apply reflecting boundary conditions for the density, velocity and 
internal energy. To represent steady conditions in the photosphere at the base of the wind, during the 
evolution of each model we apply the constraints that in the first zone above the 9 = 90° plane the radial 
velocity v r = 0, the rotational velocity remains Keplerian, and the density is fixed at p — p at all times. 
The initial conditions arc as in PSD I expect for the initial temperature profile which is described in detail 
in Section 2.3. 

To solve eqs. 1-3 we use an extended version of the ZEUS-2D code (e.g., Stone & Norman 1992). The 
equation for the internal energy, eq. 3, is solved using the operator splitting method and the backward 
Eulcr scheme. 



2.2. The radiation field and force 

The geometry and assumptions needed to compute the radiation field from the disk and central object 
are as in PSD II (sec also PSD I). The disk is flat, Keplerian, geometrically-thin and optically-thick. The 
disk photosphere is coincident with the 9 = 90° axis. We specify the radiation field of the disk by assuming 
that the temperature follows the radial profile of the optically thick accretion disk (Shakura & Sunyaev 
1973), and therefore depends on the mass accretion rate in the disk, M a , the mass of the black hole, Mbh 
and the inner edge of the disk, r* = 3rs, where r$ = 2GMbh/c 2 is the Schwarzschild radius of a black hole. 
In particular, the disk luminosity, Ld = 2r/GMBHM a /rs, where r\ is the rest mass conversion efficiency. 

The geometry of the central engine in AGNs is not well known. For simplicity, we consider the central 
engine as the most inner part of the accretion disk plus an extended corona. We refer to the corona as the 
central object. The radius of the central object is comparable with the inner radius of the disk and we 
assume that they formally are equal. 

We express the central object luminosity L* in units of the disk luminosity L* = xLn- In contrast 
to PSD I and PSD II, we allow for the situation when only some fraction of the central object luminosity 
takes part in driving a wind. We identify this fraction as the luminosity in the UV band, /uv-k*- We 
refer to the fraction of the luminosity that is responsible for ionizing the wind to a very high state as the 
luminosity in the X-ray band, /x^*- For simplicity, we assume here that this fraction of the luminosity 
does not contribute to line driving of the wind. We call the luminosity in the remaining bands, mainly 
optical and infrared, as /o,ir^*- We assume that /o,ir£* is the part of the luminosity that does not change 
the dynamics of the wind. We stress that the fraction fi is a numerical factor that we introduce here to 
parameterize the luminosity in each of those three domains. We set /opt.ir to zero in the remaining part of 
the paper as the central object is very hot. We take into account the irradiation of the disk by the central 
object, assuming that the disk re-emits all absorbed energy locally and isotropically. However we note that 
the contribution from irradiation is negligible for x ~ 1 and large radii (see eq. 6 below). 

We approximate the radiative acceleration due to lines (line force, for short) using a modified CAK 
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method. The line force at a point defined by the position vector r is 

F™*'(r) = jf M(t) ^^f d ^ (4 ) 

where / is the frequency-integrated continuum intensity in the direction defined by the unit vector n, 
and O is the solid angle subtended by the disk and central object at the point. The term in brackets is 
the electron-scattering radiation force, a e is the mass-scattering coefficient for free electrons, and M(t) is 
the force multiplier - the numerical factor which parameterizes by how much spectral lines increase the 
scattering coefficient. In the Sobolev approximation, M(t) is a function of the optical depth parameter 

\d Vl /diy (b) 

where v t h is the thermal velocity, and ^ is the velocity gradient along the line of sight, h. 

We evaluate the radiation force in four steps. First, we calculate the intensity, the velocity gradient 
in the n direction and then the optical depth parameter t. Second, we calculate the parameters of the 
force multiplier using a current value of the photoionization parameter, £ adopting results of Stevens & 
Kallman (1990). Then we calculate the radiation force exerted by radiation along h. Third, we integrate 
the radiation force over the solid angle subtended by the radiant surface. Finally, we correct the radiation 
force in the radial direction for the optical depth effects. 

Our numerical algorithm for evaluating the line force for given parameters of the force multiplier, the 
third step, is described in PSD II. For a rotating flow, there may be an azimuthal component to the line 
force even in axisymmetry. However we set this component of the line force to zero because it is always less 
than other components (e.g., PSD II). See PSD I and PSD II for further details. Below we describe our 
calculations of the force multiplier for various conditions in the wind and our treatment of the optical depth 
effects on the radiation force. 



In the disk plane at r = t\d, I(r,n) is the local isotropic disk intensity: 

Mr„) = 2e &^{i(i-(7jr)" 2 ) 



(6) 



We consider only the hot part of the disk, where the local disk temperature, To i<5 a few 10 3 K. We then 
assume that all disk photons can contribute to the line force. The second term in the curly brackets 
corresponds to the contribution from the irradiation of a disk by a central object. For large radii and x ~ 1, 
the contribution from the irradiation is negligible compared to the intrinsic disk intensity (the first term in 
the curly brackets). 

For radiation from the central object, the intensity may be written as: 

= ™- (7) 



We account for the fact that the central radiation consists of the two distinct spectral components by using 
the parameter where i — X or UV. We assume that the optical depth effects in the 9 direction are 
negligible for the central radiation in the two bands. Note that the 9 component of the radiation force due 
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to the central object is negligible in our computational domain because the domain is far from the object 
(r* CrJ. In the radial direction however, the column density can be high and the optical depth effects can 
be significant. We estimate the optical depth, t\ between the central source and a point in a wind from: 



n = ( Kipdr, (8) 
Jo 



where Ki is the absorption coefficient representative for the i band, and r is the distance from the central 
source. We can treat the central object as a point source because our calculations are for r*/r <C 1. Then 
we calculate the radial component of the central radiation force due to lines from 

(Fl ad \ (r) = /uvexp(-r uv ) / M(i) ( n /Mr,n)dQ \ ^ (g) 



n, \ c 

where 17* is the central object solid angle. We also take into account the optical depth effects on the radial 
component of the electron-scattering force: 

{F rad,e )r (r) = {/^ e xp(-T UV ) + / X CX P (-T X )} I L /M^n)dQ \ _ (1Q) 



The attenuation of the X-ray radiation is calculated using k x = 40 g 1 cm 2 for £ < 10 5 and k x = 0.4 g 1 cm 2 
for £ > 10 5 , while the attenuation of the UV radiation is calculated using k uv = 0-4 g _1 cm 2 for all £. 

A disk wind may also be optically thick to the UV radiation emitted by the disk. In particular, the 
flux along the radial direction near the disk can be significantly reduced due to high column density. We 
approximate this effect by multiplying the radial component of the disk radiation force by the attenuation 
factor exp(— ruv)- 

To calculate the force multiplier, we adopt the CAK analytical expression modified by Owocki, Castor 
& Rybicki (1988, see also PSD I) 



M(t) = kt~ 



(l+T V 1 -") - 1 
V- 1 - ~ ' max J - 1 - 



T (l"«) 
' max 



(ii) 



where k is proportional to the total number of lines, a is the ratio of optically-thick to optically-thin lines, 
Tmax = tr/max and i] max is a parameter determining the maximum value, M max achieved for the force 
multiplier. Equation 11 shows the following limiting behavior: 

lim M(t) = kr a (12) 
lim M(t) = M max , (13) 

where M max = k(l — a)rim ax . The maximum value of the force multiplier is a function of physical 
parameters of the wind and radiation field. Several studies have showed that M max is roughly a few 
thousand for gas ionized by a weak or moderate radiation field (e.g., CAK; Abbott 1982; Stevens & Kallman 
1990; Gayley 1995). As the radiation field becomes stronger and the gas becomes more ionized the force 
multiplier decreases asymptotically to zero. 

The line force, in particular the parameters of the force multiplier depend on the ionization of the wind 
and the spectral energy distribution (SED) of the radiation field. Self-consistent calculations of the line-force 
for given wind conditions and the radiation require detailed calculations of the wind photoionization 
structure (e.g., CAK; Abbott 1982; Puis et al. 2000). The total line force includes contributions from > 10 5 
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lines from many species. Such calculations in connection with the 2-D, time-dependent hydrodynamic 
calculations are not feasible. We therefore start by adopting the analytical formulae for the force multiplier 
due to Stevens & Kallman (1991). They studied the effects of X-ray ionization on the the radiative force 
experienced by the stellar wind in a massive X-ray binary (MXRB). Wind conditions in MXRBs differ 
somewhat from those in QSOs. In particular, the X-ray radiation in MXRBs can be well represented by a 
10 KeV bremsstrahlung spectrum while the spectrum in QSOs is better fitted by power laws with different 
spectral indices at different spectral bands (e.g., Zheng et al. 1997; Laor et al. 1997). Nevertheless we adopt 
the Stevens & Kallman results as a first order approximation because the Compton temperature and the 
ionization structure are similar (Kallman & McCray 1982). 

Using the photoionization code XSTAR, Stevens & Kallman found that the line force due to the 
radiation from the primary can be parameterized in terms of the CAK force multiplier. Their results show 
that the line force decreases sharply with increasing photoionization parameter: 

*=^, (14) 
n 

where Tx is the local X-ray flux, n is the number density of the gas (= p/(m p fi) , where m p is the proton 
mass, and [i is the mean molecular weight). Stevens & Kallman also found some simple analytical fits to 
their results that allowed them to express k and i] max as functions of £: 

k = 0.03 + 0.385 cxp(-1.4£ - 6 ), (15) 



and 

r 6.9 exp(0.16 £ ' 4 ) for log 10 £ < 0.5 

logio Vmax = I (16) 

[ 9.1cxp(-7.96 x IQ~ 3 for log 10 £ > 0.5 
The parameter a = 0.6 and does not change with £. 

These expressions for the parameters of the force multiplier predict that M max increases gradually 
from ~ 2000 to 5000 as £ increases from to ~ 3 and then drops to ~ 0.1 at £ = 1000. The line force 
becomes negligible for £ ^ 100 because then M max ~ 1. 

To calculate the photoionization parameter, we first need to know the local X-ray flux Tx (see eq. 14). 
We estimate Tx assuming that the source of all X-rays is a point-like central object with the luminosity in 
the X-ray band, Lx = fxL*. Making the above simplifications, we can express the local X-ray flux as 

Tx = cM-rx)^. (17) 



2.3. Radiation heating and cooling 



We calculate the gas temperature assuming that the gas is optically thin to its own cooling radiation 
and that the abundances are cosmic (Withbroe 1971). The net cooling rate depends on the density, /?, 
the temperature, T, the ionization parameter £, and the characteristic temperature of the X-ray radiation 
Tx- It is then possible to fit analytical formulae to the heating and cooling rate obtained from detailed 
photoionization calculations for various Tx, T, and £. For example, Blondin (1994) found that for a 10 keV 



- 9 - 



bremsstrahlung spectrum, his fits typical agree with a detailed calculation (Blondin et al. 1990) to within 
25%. Using Blondin's result we can express the net cooling rate C in equation 3 by 

pC = n 2 (G C o mp ton + G x - Lb,i) erg cm 3 s -1 , (18) 

where Gcompton is the rate of Compton heating/cooling, 

Gcompton = 8.9 x 10~ 36 £ (T x - T) erg cm" 3 s" 1 , (19) 

Gx is the net rate of X-ray photoionization heating-recombination cooling 

Gx = 1.5 x 1(T 21 e 1/4 T- X I 2 {1 -T/T x ) erg cm" 3 s" 1 (20) 

and L^i is the rate of bremsstrahlung and line cooling 

L w = 3.3 x 10~ 27 T 1 / 2 + [1.7 x 10- 18 exp(-1.3 x 10 5 /T^T' 1 ' 2 + 10" 24 ]<5 erg cm" 3 s" 1 . (21) 

The parameter 8 is introduced to control line cooling, 8 — 1 represents optically thin cooling and 8 < 1 
represents the case when lines become optically thick and cooling is reduced. 

As we mentioned in Section 2.1, we use the same initial and boundary conditions as in PSD I. However 
we here also calculate the evolution of the internal energy. We specify the boundary conditions and the 
initial conditions for e as follows. For the temperature at the base of the wind, in the first grid zone above 
the 9 = 90° plane, we adopt the solution for the steady state disk. Adopting the assumption that the local 
disk intensity is a black body and using eq. 6, the disk temperature is: 

T D {r D ) -- 

Then we calculate the internal energy from 

e a {r D ) 

where p is the density in the first grid zone above the 9 = 90° plane. We consider models with \i — 1. We 
initialize the internal energy at other locations by assuming that the gas temperature at a given radius is 
the same as at the disk temperature (i.e., T(r, 9) = To(r sin (9) for all 9) and using initial value of the density 
(see PSD I). Consistent with the boundary conditions for other gas quantities (see Sect. 2.1) we apply the 
constraint that in the first zone above the 9 = 90° plane the internal energy e = e at all times. We find 
only the initial, transient evolution of the wind is affected by this choice of initial temperature conditions. 

Our treatment of radiative heating and cooling is most likely valid in a low density and high 
photoionization parameter regime. For a high density and low photoionization parameter regime our 
treatment will likely underestimate the temperature because we oversimplify line cooling by not properly 
taking into account optical depth effects. In dense regions, our scheme therefore may yield a very low gas 
temperature (e.g., T <C 1000 K). In such cases we set the lower limit for the temperature assuming that 
the gas is in local thermodynamical equilibrium, i.e., the gas temperature equals the local disk radiation 
temperature. Specifically, if at the (r, 9) wind point the temperature from eq. 3 is less than Tr>(rsin9) then 
we replace it with Toirsinff). 



(7rI D (r D )/a)^ 4 . (22) 

= PokT{r D ) 
pm p (~f~- 1)' 



3. Results 
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We specify our models by several parameters. In all our calculations we assume the mass of the 
non-rotating black hole, Mbh = 10 8 M and the rest mass conversion efficiency, 77 = 0.06 typical for AGN. 
To determine the radiation field from the disk, we assume the mass accretion rate M a = 1.8 M Q yr _1 . These 
system parameters yield the disk Eddington number, T D = 0.5, the disk inner radius, = 8.8 x 10 13 cm 

and the orbital period at the inner edge of the disk, r = \J G m bh = 7.22 x 10 3 sec. We set the offset of 
the computational domain, z = Hjj w 3Td that corresponds to the half width of the radiation pressure 
dominated disk at large radii (e.g., Shakura & Sunyaev 1973). The radiation field from the central engine 
is specified by the additional parameters: x — 1, /uv = 0.5, fx = 0.5. As we discussed above the SED 
of the ionizing radiation is not well known, our choice of values for /uv and fx is guided by the results 
from Zheng et al. (1997) and Laor et al. (1997). To calculate the line force we adopt the force multiplier 
parameter a = 0.6. We calculate the remaining parameters of the multiplier, i.e., k and rj m ax as functions 
of the photoionization parameter, £ using cqs 15 and 16. The force multiplier depends only formally on the 
thermal speed, vth which we set to 20 km s _1 , i.e., the thermal speed of a hydrogen atom at the temperature 
of 25000 K (Stevens & Kallman 1990). To calculate the gas temperature, we assume the temperature of the 
X-ray radiation, Tx = 10 KeV and the line cooling parameter 5 = 1. 

We adopted a generalized CAK method to calculate a disk wind driven by the line force (PSD II). The 
CAK method has been developed for OB stars with photospheres dominated by gas pressure (the inner 
atmosphere is in hydrostatic equilibrium), and they radiate most of their energy in the UV band. We chose 
the radial extent of the computational domain (i.e., 100 r* < r < 1000 r*) keeping in mind that the disk 
atmosphere should satisfy these two conditions, however some compromises had to be done. For example, 
the disk temperature at 100 r* = 8.8 x 10 15 cm is 8200 K while at 1000 r* = 8.8 x 10 16 cm is 1500 K using 
our disk parameters. We realize that according to the a-disk model the disk emits mostly in the UV band 
down to the radius of ~ 10 r* where its effective temperature is 50000 K. However the disk at this small 
radius is likely radiation dominated (e.g., Shakura & Sunyaev 1973; Svensson & Zdziarski 1994). On the 
other hand, the disk atmosphere at radius 1000 r* and the effective temperature of 1500 K, radiates most of 
its energy in the infrared band. Nevertheless we chose a large value for the outer edge of the computational 
domain to make sure that the domain includes most of the acceleration zone of the disk wind launched from 
the inner disk and the wind velocity at the outer edge well represents the wind terminal velocity. 

Figure 2 shows the instantaneous density, temperature and photoionization parameter distributions 
and the poloidal velocity field of the model. After 12.6 years (~ 5.5 x 10 4 r), the disk material fills the 
grid for 9 ^ 70° and remains in that region for the rest of the run, i.e., over next 47.7 years. Figure 2 shows 
results after 14.6 years. Although the flow is time-dependent the gross properties of the flow (e.g., the mass 
loss rate and the radial velocity at the outer boundary), settle down to steady time- averages over timescales 
on the order of 3 years. Our calculation follows (i) a hot, low density flow in the polar region (ii) a dense, 
warm and fast equatorial outflow from the disk, (iii) a transitional zone in which the disk outflow is hot and 
struggles to escape the system. 

In the polar region, the density is very small and close to the lower limit that we set on the grid, i.e., 
Pmin = 10~ 20 gcm~ 3 . The line force is negligible because the matter is highly ionized as indicated by a very 
large photoionization parameter (~ 10 7 ). The gas temperature is close to the temperature of the X-ray 
radiation, again indicative of highly ionized gas. The matter in the polar region is pulled by the gravity 
from the outer boundary and it is an artifact of the boundary conditions. Overall this region of the very 
low density is not relevant to our analysis as it has no effect on the much denser disk flow. 

The equatorial region is distinctly different. In the inner part of the disk (i.e., for r ~ r^), the density 
at the wind base is high, ~ 10~ 13 g cm~ 3 . Thus the photoionization parameter is low despite the strong 
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central radiation. However as the flow from the inner part of the disk is accelerated by the line force its 
density decreases and the gas temperature and the photoionization parameter increase. Subsequently the 
gas becomes fully ionized and loses all of driving lines before it reaches the escape velocity and therefore 
falls on the central object/inner disk. Although this gas does not produce a wind, its primary effect is to 
shield the gas at larger radii. Thus the wind consists of gas accelerated by the line force at larger radii, in 
fact driving of the disk wind extends over all radii at which the intrinsic disk radiation is large enough to 
launch gas. 

The poloidal velocity (Fig. 2c) shows that the gas streamlines arc perpendicular to the disk over some 
height that increases with radius. The streamlines then bend away from the central object and converge. 
The region where the flow is moving almost radially outward is associated with a high- velocity, high density 
stream. This fast stream contributes ~ 100% to the total mass loss rate, Mw = 0.5 M yr -1 . We note that 
the mass loss rate can increase by a factor of a few when a knot is crossing the outer boundary. 

The fast stream is variable. In the upper envelope of the disk wind there is large velocity shear 
between the higher density fast stream moving outward and the lower density fast gas moving inward. Our 
simulations show that this shear gives rise to Kelvin-Hclmholtz instabilities. The instabilities generate 
knots that propagate along the fast stream at ~ 10000 km s _1 . Figure 2a shows an example of such a knot 
at r ~ 750 r* and z ~ 200 r„, in the figure coordinates. The density contrast between the knot and the fast 
stream is ~ 2 orders of magnitude. The knots generator is episodic and is inherent to the fast stream. We 
observed the generation of 19 knots over 60.4 years. In other words, a knot is produced every ~ 3 years. 

Figure 3 presents a sequence of density maps showing time evolution of the outflow from Figure 2 
after 13.3, 14.6 and 16.47 years, left, middle and righ panel respectively. Note that the middle panel 
from Figure 3 is the same as the top left panel from Figure 2. Figure 3 well illustrates variability in the 
outflow, in particular, generation of a knot in the fast stream at r ~ 450 and z ~ 100 (left panel), the 
well- formed knot at r ~ 750 r* and z ~ 200 r* (middle panel), and the time when the knot left the gird 
at r <~ 950 r* and z <~ 300 r*, and the fast stream is fairly smooth in between episodes with knots (right 
panel) . The time dependence in the region close to the disk, bounded by the fast stream resembles outward 
propagation of a wave at the begin of evolution but becomes more complex with time as elements of the 
flow move both upwards and downwards. 

Figure 4 presents the run of the density, radial velocity, photoionization parameter and column density 
as a function of the polar angle, 9 at the outer boundary, r Q — 8.8 x 10 16 cm from Figure 2. The column 
density is given by: 



The gas density is a very strong function of angle for 9 between 90° and 65°. The density drops by 8 
orders of magnitude between 9 = 90° and 9 ~ 89°, as might be expected of a density profile determined by 
hydrostatic equilibrium. For 70° Ss 9 89°, the wind domain, p varies between 10~ 19 and 10~ 17 g cm~ 3 . 
For 9 Ss 70°, density decreases gradually to so low a value that it becomes necessary to replace it by 
the numerical lower limit p m in- The radial velocity at 1000 r* = 8.8 x 10 16 cm has a broad peak for 
72° <0< 82° with the maximum of 15000 km s _1 at 9 ~ 75. On both side of the peak the velocity is close 
to zero. 

The photoionization parameter is very high, ~ 10 7 , for 9 ^ 70° because of the very low density. 
However it drops by 15 orders of magnitude between 70° Ss 9 75° and stays at a very low level for 
9 > 75°. The column density changes less dramatically with angle. Between the pole and 9 ^ 67°, N H is 




(24) 
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less than 10 22 cm 2 . For > 67°, Nh increases gradually with 8, it reaches value of 10 24 cm 2 at 9 ~ 80°. 
For 80° < (9 < 89°, N H varies between 10 24 cm" 2 and 10 26 cm" 2 . 

Our model shows that (i) the intrinsic disk radiation can launch a wind and (ii) the central object 
radiation can accelerate the disk wind to very high velocities. We have checked how these main results 
are sensitive to the model parameters. For a fixed disk atmosphere and central radiation source, the most 
important parameter of our model is the optical depth for the X-rays from the central object. Therefore we 
have run several tests with various X-ray opacities for £ < 10 5 : kx = 0.4 g _1 cm 2 and kx = 4 g _1 cm 2 . 
The former case is most conservative because it corresponds to the case when the X-ray optical depth is 
reduced to the Thomson optical depth. We found that in both test runs the X-ray radiation is significantly 
attenuated so the line force can launch a disk outflow. However the flow velocity never exceeds the escape 
velocity because the X-rays fully ionize the gas close to the disk and produces a hot corona with complex 
velocity field but which does not escape the system. We conclude then that the disk atmosphere can 'shield' 
itself at least to the extent that the local disk radiation can launch gas off the disk photosphere. 

The two test runs also illustrate how robust is our second result: the central radiation accelerates the 
wind to very high velocities. For this to happen the column density, Nn must be large enough to reduce 
the X-ray radiation but too not large to reduce the UV flux in the radial direction. In other words, this 
requires r x > 1 and at the same time, tjjv *C 1. 

The test run with kx — 0.4 g _1 cm 2 corresponds to the situation where rx = t\jv- Thus the disk gas 
that is sufficiently shielded from the X-rays is also shielded from the UV photons from the central object. 
Such gas can only be lifted from the disk surface by the disk UV radiation but fails to gain momentum in 
the radial direction. 

The fact that the disk wind can be launched without external shielding material implies that our 
solution can depend on the location of the inner computational radius, rj. We have run several tests with 
Ti = 50r» and n = 200r» to check this (the remaining model parameters were as in the models shown 
in Figure 1 and 2). We found that the location of the inner edge, affects the properties of the wind 
but not the fact the wind is produced. For example, the wind opening angle, u>, is 25°, 15° and 8° for 
ri = 50, r*, lOOr* and 200?%, respectively. The maximum radial velocity at the outer boundary decreases 
from 20000 km s _1 through, 15000 km s _1 to 5000 km s _1 for above sequence of decreasing r». The mass 
loss rate also decreases with n from 0.6 M Q yr -1 through 0.5 M Q yr _1 to 0.2 M Q yr _1 . 



4. Discussion 

Our calculations show that the UV emitting accretion disk can launch a wind that will shield itself from 
the strong ionizing radiation emitted from the central object. The self-shielding is quite robust because the 
external radiation does not penetrate the disk down to the region where the line force becomes important 
and pushes matter away from the photosphere. We illustrated the robustness of the self-shielding by 
changing relevant model parameters. We also note that a similar effect was observed by Stevens & Kallman 
(1990) and Stevens (1991). They studied the effects of X-ray ionization on the radiative force experienced 
by the stellar wind in a MXRB. They calculated the detailed photoionization structure of the line-driven 
wind for various X-ray luminosities without optical depth effects (Stevens & Kallman 1990) and with optical 
depth effects (Stevens 1991). Additionally, they considered irradiation at various directions including the 
normal to the mass losing photosphere. In none of the models the irradiation penetrated below the so-called 
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wind critical point where, in the CAK type models, the wind mass loss rate is determined. They found, as 
we did, that the irradiation can severely decrease the wind velocity in the supersonic portion of the flow. 

We model here a line-driven wind from a disk that is flat, Keplcrian, geometrically-thin and 
optically-thick. For the temperature at the base of the wind, in the first grid zone above the 9 = 90° plane, 
we adopt the solution for the steady state a disk. We calculate the local disk intensity as if it emits as a 
black body. In all our calculations, we treat the regions for 9 ~ 90° as if they are gas pressure dominated. 
We find that our results depend on the inner radius of the computational domain. This dependence is 
likely an artifact of one or more of our assumptions. For example, for very high luminosities the inner part 
of the disk is dominated by the radiation pressure. For the parameters adopted here, namely V — 0.5, 
the inner radius of a gas pressure dominated disk is ~ 10 16 cm (e.g., Svensson & Zdziarski 1994). The 
structure of the radiation dominated disks is not well known. We expect, however, that incorporating the 
radiation dominated part of the disk into our calculation will change the condition of the base of the wind, 
for example, the gas density and the UV flux might be reduced. These changes might affect the solution for 
the wind, in particular, the wind mass flux. We plan to explore these issues in a future paper. 

To calculate the disk surface temperature and intensity, we took into account the irradiation of 
the disk by the central object. However the contribution from the irradiation is negligible compared to 
the disk intrinsic temperature because we considered a flat disk at large radii. Additionally, the high 
column density of the disk wind, that we found in our models, implies that the irradiating flux will be 
significantly attenuated by the wind before it reaches the disk surface. We anticipate that the last effect 
will significantly reduce disk irradiation regardless of the shape of the disk surface - flat or flaring. Detailed 
NLTE photoionization calculations are required to determine what fraction of the central radiation will 
reach the surface of a mass losing disk. A complete treatment of disk irradiation is difficult also because it 
is not obvious a priori whether the radiation incident on the disk is completely thermalized and re-radiated 
isotropically or it is scattered off the disk atmosphere or both. 

Our calculations show that a line-driven disk wind model offers a promise of explaining outflows in 
AGNs. However they also illustrate some problems with this model. For example, as dKB pointed out 
in discussing MCGV's model, a very small radius at which the disk wind is launched also implies a very 
small size to BELR because BELR lies inside or is cospatial with BALR. The small scales in turn imply 
short crossing time of the BALR, of order of a month or so, and it is difficult to understand that the highly 
complex kinematical structures that BALs often exhibit do not appear to vary on timescales of 10 years 
(Barlow 1994). 

Our axisymmetric 2D models show that the disk wind is unsteady and generates dense knots every 
3 years. These knots will correspond to rings in 3D. However in fully 3D calculations the rings may break 
down to spirals or clouds or both. The breaking down of the rings will change the density contrast between 
them and the rest of the wind. Thus it is not clear if we would be able to see any spectral signature 
of density fluctuations in the wind. Detail calculations of line profiles are required to check this point. 
Additionally, detailed photoionization calculations are required to check if the full range of ions observed to 
show BAL profiles can be explained: BALs from ions with ionization potentials as low as of O HI or lower 
and as high as of O VI. 

To produce a fast wind the ratio between rx to tuv is very important. The low ratio gives a slow disk 
wind whereas the high ratio gives a fast wind. This result is consistent with the observed anti-correlation 
for QSOs between the relative strength of the soft X-ray flux and the CIV absorption equivalent width 
(e.g., Brandt, Laor & Wills 2000). However there is an upper limit for the X-ray attenuation, or the column 
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dcnsity between the X-ray source and wind, namely for Nu ^ 10 24 cm -2 (jx >> 1) the gas is well shielded 
from the X-rays but at the same time it is also shielded from any other radiation from the central engine 
including the UV radiation. In this case, the line force from the central engine could be so much reduced 
that it may not accelerate the gas to high radial velocities and no strong wind will be produced. 

Some of the emission-line properties of QSOs could be explained as a correlation between luminosity 
and the slope of the ionizing spectrum, i.e., lower luminosity objects have harder spectra (e.g., Boroson & 
Green 1992). Using our parameterization of the AGN radiation, this QSO's property can be represented by 
increasing fx in expense of /uv when we reduce the luminosity. Such a choice of model parameters will 
increase the strength of X-ray ionizing radiation in comparison to the UV driving radiation. Consequently 
a disk wind should be weaker and slower as some of our test runs indicate. 

Boroson & Green (1992) argue that the dominant source of variation in the observed properties of low 
redshift QSOs is not driven by external orientation but rather by the fraction of the Eddington luminosity 
at which the object is emitting (our parameter To). We plan to examine the parameter space of our models 
to define the major trends in disk wind behavior that will help us to explain the observational trends found 
in various AGNs. 



5. Conclusions 

We have studied radiation driven winds from luminous accretion disks using numerical methods to 
solve the two-dimensional, time-dependent equations of hydrodynamics. In so doing we have accounted for 
the radiation force due to spectral lines using a generalized multidimensional formulation of the Sobolev 
approximation. Additionally we have taken into account the effects of the strong central radiation on the 
wind photoionization structure and thermodynamics. 

We find that the local disk radiation can launch a wind from the disk despite strong ionizing radiation 
from the central object. The central radiation may overionize the supersonic portion of the flow and severely 
reduce the wind velocity. To produce a fast disk wind the wind X-ray opacity must be higher than the UV 
opacity by <; 2 orders of magnitude. 

Our calculations of a wind from a disk accreting onto a 10 s M Q black hole at the rate of 1.8 M© yr" 1 
show that the radiation force from an accretion disk can launch a self-shielding wind from a radius of 
55 10 16 cm while the strong UV radiation from the central object can radially accelerate the disk wind to 
velocities ~ 15000 km s _1 , for the X-ray opacity of 40 g" 1 cm 2 . The disk wind domain is intrinsically 
unsteady and its covering factor is ~ 0.2. The wind mass loss is 0.5 M Q yr" 1 which is a significant fraction 
of the mass accretion rate. The strong X-ray radiation from a central object completely ionizes the polar 
region and only a thin layer above the upper envelope of the disk wind. The disk wind immediately below 
the upper envelope can be characterized as a fast, high-density stream which is reminiscent of the PSD 
disk wind solution dominated by the driving radiation from a bright central object. The column density of 
the fast stream is between 10 22 cm" 2 and 10 24 cm" 2 so the stream is optically thin to the UV radiation 
and this is precisely why it can be accelerated to high velocities. On the other hand, the part of the wind 
bounded by the fast stream, closer to the disk, is slow and is reminiscent of the PSD solution dominated by 
the driving radiation from the disk. The column density of the slow wind is > 10 24 cm" 2 so it is completely 
shielded from the central radiation, both the X-rays and the UV photons. 
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Fig. 1. — The framework of our two-dimensional hydrodynamical calculations for a line-driven disk wind. 
The drawing is not to scale. We assume the disk is fiat, Keplerian, geometrically thin, and optically thick. 
The radiation pressure dominated disk is represented by the blue regions while the gas pressure dominated 
disk is represented by the red regions. The black hole is represented by the black circle while the X-ray 
source, is marked by the yellow region. The dashed line perpendicular to the disk is the disk rotational axis. 
Our computational domain is marked by solid lines. The = 90° axis, the dashed line perpendicular to the 
rotational axis, is located above the disk midplane, the dot-dashed line. The offset of the 8 = 90° axis from 
the disk midplane is given by the disk pressure scale height. See Section 2 for more details. 
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Fig. 2. — The top left panel is a color density map of the AGN disk wind model, described in the text. 
The top right panel is a color gas temperature map of the model while the bottom left panel is a color 
photoionization parameter map. Finally, the bottom right panel is a map of the velocity field (the poloidal 
component only). In all panels the rotation axis of the disk is along the left hand vertical frame, while the 
midplane of the disk is along the lower horizontal frame. 
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Fig. 3. — A sequence of density maps for the model from Figure 2 after 13.3, 14.6, and 16.47 years (left, 
middle, and right panel). Note both the time-dependent fine structure near the base of the wind and the 
time-dependent large scale structure associated with the fast stream at the polar angle ~ 75° 
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Fig. 4. — Quantities at the outer boundary in our model (the top left panel of Figure 2 and the middle panel 
of Figure 3). The ordinate on the left hand side of each panel refers to the solid line, while the ordinate on 
the right hand side refers to the dotted line. The column density, Nh is calculated along the radial direction. 
The radial velocity peaks at 9 « 74° while the density at 9 w 71°. The column density is less than 10 24 cm~ 2 
for 9 < 80°. 



